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Abstract 

In striated muscle, the number of capillaries containing moving red blood 
cells increases with increasing metabolic demand. This phenomenon, termed 
capillary recruitment, has long been recognized, but its mechanism has been 
unclear. Here, a theoretical model for metabolic blood flow regulation in a 
heterogeneous network is used to test the hypothesis that capillary recruitment 
occurs as a result of active control of arteriolar diameters, combined with 
unequal partition of hematocrit at diverging microvascular bifurcations. The 
network structure is derived from published observations of hamster cremas- 
ter muscle in control and dilated states. The model for modulation of arterio- 
lar diameters includes length-tension characteristics of vascular smooth muscle 
and responses of smooth muscle tone to myogenic, shear-dependent, and met- 
abolic stimuli. Blood flow is simulated including nonuniform hematocrit dis- 
tribution. Convective and diffusive oxygen transport in the network is 
simulated. Oxygen-dependent metabolic signals are assumed to be conducted 
upstream from distal vessels to arterioles. With increasing oxygen demand, 
arterioles dilate, blood flow increases, and the numbers of flowing arterioles 
and capillaries, as defined by red blood cell flux above a small threshold value, 
increase. Unequal hematocrit partition at diverging bifurcations contributes to 
recruitment and enhances tissue oxygenation. The results imply that capillary 
recruitment, as observed in the hamster cremaster preparations, can occur as 
a consequence of local control of arteriolar tone and the resulting nonuniform 
changes in red blood cell fluxes, and provide an explanation for observations 
of sequential recruitment of individual capillaries in response to modulation 
of terminal arteriolar diameter. 



Introduction 

Proper distribution of blood flow within organs is essen- 
tial for the matching of oxygen and nutrient supply to tis- 
sue demand. In the microcirculation, local regulation of 
blood flow is achieved by contraction and relaxation of 
the vascular smooth muscle (VSM) surrounding the arte- 
rioles and small arteries. Changes in VSM tone are elic- 
ited in response to changes in several factors, including 
levels of oxygen and other metabolites, intraluminal pres- 
sure, and luminal wall shear stress (Duling et al. 1987). 
These vasoactive stimuli not only act locally but also 
induce responses that are propagated upstream, causing 



vessels feeding the site of the stimulus to constrict or 
dilate. The resulting coordinated control of VSM tone in 
vessels proximal to the affected site contributes to blood 
flow regulation in response to changes in the metabolic 
demand of the tissue (Segal et al. 1989; Berg et al. 1997; 
Cohen et al. 2000). 

Blood flow rate is regulated over a particularly wide 
range in skeletal muscle, which experiences large varia- 
tions in metabolic demand (Saltin et al. 1998; Segal 
2005). How perfusion is matched to demand in skeletal 
muscle has not been fully resolved. Starting with the work 
of Krogh (1919a,b), it has been proposed that one of the 
contributing mechanisms is capillary recruitment, 

2013 | Vol. 1 | Iss. 3 | e00050 
Page 1 



© 2013 The Authors. Physiological Reports published by Wiley Periodicals, Inc. on behalf of 
the American Physiological Society and The Physiological Society. 

This is an open access article under the terms of the Creative Commons Attribution License, 
which permits use, distribution and reproduction in any medium, provided the original work is properly cited. 



Metabolic Flow Regulation and Capillary Recruitment 



B. C. Fry. ef al. 



whereby vessels that are not flowing at rest commence 
flowing in response to increased metabolic demand 
(Klitzman et al. 1982; Berg et al. 1997). The mechanisms 
of recruitment, however, remain incompletely understood, 
and even the occurrence of recruitment is subject to 
debate (Clark et al. 2008). 

Early investigators considered that capillary recruitment 
resulted from control of flow in individual capillaries by 
"pre-capillary sphincters" (Krogh 1919a,b). This concept 
has persisted in textbooks. However, precapillary sphinc- 
ters are not observed in skeletal muscle (Gorczynski et al. 
1978; Lindbom and Arfors 1985; Johnson 1995; Segal 
2005). Instead, active control of flow in capillaries of skel- 
etal muscle has been observed to occur only in vessels 
down to the level of terminal arterioles (Delashaw and 
Duling 1988; Berg 1995). These observations led to the 
concept that the smallest unit for control of blood flow 
in skeletal muscle is the microvascular unit, consisting of 
a terminal arteriole and group of capillaries fed by it 
(numbering about 20) (Delashaw and Duling 1988; Emerson 
and Segal 1997; Lo et al. 2003). 

According to this concept, recruitment would be 
expected to occur at the level of microvascular units, not 
individual capillaries. However, Lindbom and Arfors 
(1985) observed a sequential increase in the number of 
RBC-perfused capillaries fed by a single terminal arteriole 
during vasodilation elicited by a decrease in ambient oxy- 
gen levels. Their observations imply that recruitment of 
individual capillaries can occur in response to dilation of 
terminal arterioles. Lindbom and Arfors (1985) attributed 
this behavior to the existence of a yield shear stress in 
blood. However, the yield shear stress refers to bulk blood 
flow and does not apply to red blood cells (RBCs) flow- 
ing in narrow capillaries. An alternative potential explana- 
tion is that nonlinear rheological effects associated with 
RBC motion in nonuniform capillaries cause flow cessa- 
tion. A previous quantitative analysis (Secomb and Hsu 
1995) showed that such behavior would occur only at 
flow velocities much lower than those seen in experiments 
(Lindbom and Arfors 1985), and cannot account for the 
observed recruitment. Thus, the mechanism of recruit- 
ment at the level of individual capillaries in response to 
modulation of arteriolar diameters has remained unclear. 

The regulation of blood flow in the microcirculation 
involves a number of interacting physical and biological 
processes, occurring in a network with a complex and 
heterogeneous structure. Due to technical limitations, it is 
not possible to measure many of the relevant variables in 
individual microvessels. Therefore, theoretical models 
have played an important role in the investigation of mi- 
crocirculatory phenomena (Secomb et al. 2008). Such 
models provide a framework for integrating available 
information and quantitatively testing hypotheses regard- 



ing microcirculatory function. The effects of varying the 
number and location of perfused capillaries on tissue 
oxygenation in skeletal muscle have been studied using 
oxygen transport simulations (Lo et al. 2003; Goldman 
et al. 2004). These models did not address the mecha- 
nisms determining the distributions of flowing vessels. 
Previous theoretical models of flow regulation (Ursino 
et al. 1998; Cornelissen et al. 2002; Arciero et al. 2008; 
Carlson et al. 2008; Kleinstreuer et al. 2008) have gener- 
ally considered vascular networks as sets of compartments 
connected in series, where each compartment contains 
multiple identical vessels of a given type, connected in 
parallel. Such models do not lend themselves to studying 
recruitment, which necessarily involves differences in 
behavior among capillaries, with some ceasing to flow 
while others continue flowing. 

In microvascular networks with heterogeneous struc- 
tures, experiments have demonstrated unequal partition 
of discharge hematocrit (H D ) at diverging bifurcations, 
such that the daughter vessel with the larger flow rate 
tends to receive a larger hematocrit than the other daugh- 
ter vessel (Schmid-Schonbein et al. 1980). This phenome- 
non, also referred to as phase separation, creates a 
nonuniform distribution of hematocrit throughout the 
network. A recent study (Roy et al. 2012) examined meta- 
bolic flow regulation in a heterogeneous network derived 
from the rat mesentery, assuming that the metabolic sig- 
nal was derived from RBCs. Including phase separation in 
the model led to different behavior than was predicted if 
uniform hematocrit was assumed. This shows the impor- 
tance of considering effects of phase separation when ana- 
lyzing flow regulation in heterogeneous networks. 

The objective of this study is to test the hypothesis 
that capillary recruitment occurs in heterogeneous micro- 
vascular networks as a consequence of local blood flow 
regulation by changes in arteriolar VSM tone, in combi- 
nation with unequal hematocrit partition at diverging 
microvessel bifurcations. This concept was proposed 
previously (Honig et al. 1980; Schmid-Schonbein and 
Murfee 2008), but it was not demonstrated that this 
mechanism provides a quantitative explanation for the 
observed capillary recruitment. Here, theoretical models 
are used to simulate blood flow, oxygen transport, and 
flow regulation in a heterogeneous microvascular network 
structure that is based on previously published experi- 
mental results in the hamster cremaster muscle (Berg 
1995). The simulation of blood flow includes effects of 
phase separation. Oxygen transport is simulated using a 
method that takes into account the effects of all sur- 
rounding vessels on the oxygen level at each point in the 
tissue (Secomb et al. 2004). The theoretical model for 
metabolic flow regulation is based on modulation of 
arteriolar diameters according to the length-tension char- 



2013 | Vol. 1 | Iss. 3 | e00050 
Page 2 



© 2013 The Authors. Physiological Reports published by Wiley Periodicals, Inc. on behalf of 
the American Physiological Society and The Physiological Society. 



B. C. Fry. ef al. 



Metabolic Flow Regulation and Capillary Recruitment 



acteristics of VSM (Carlson and Secomb 2005; Arciero 
et al. 2008; Carlson et al. 2008; Roy et al. 2012). 
Responses of VSM tone to myogenic, shear-dependent, 
and metabolic stimuli are included. 

Methods 

Network structure 

Most previous models of blood flow regulation have 
described the microcirculation in terms of a set of repre- 
sentative vessel segments (Ursino et al. 1998; Cornelissen 
et al. 2002; Arciero et al. 2008; Carlson et al. 2008). In 
such models, the vasculature is divided into several classes 
of vessels, for example, artery, large arteriole, small arteri- 
ole, capillary, small venule, large venule, and vein. The 
vessels in a particular class are assumed to be equivalent 
and arranged in parallel, with the different classes 
arranged in series. This assumption simplifies computa- 
tions, but does not take into account structural heteroge- 
neities present in actual microvasculature. In the study of 
Roy et al. (2012), the heterogeneity of network structure 
was represented by replacing the small arterioles, capillar- 
ies, and small venules in the representative segment 
model with a realistic network structure derived from 
experimental observations in the rat mesentery. A similar 
approach is used in this study, as shown in Figure 1. 
Here, the network structure is derived from observations 
of the hamster cremaster muscle (Berg 1995; Berg et al. 
1997), in which the microvasculature was imaged and 
mapped in a control state and in a pharmacologically 
induced maximally dilated state, and RBC fluxes were 
measured in the arterioles and capillaries. Some vessels 



had no observable RBC flux in the control state, but 
observable RBC flux in the maximally dilated state, indi- 
cating the occurrence of recruitment. All the vessels with 
observable RBC fluxes in the control state or the dilated 
state or both were included in the network maps for the 
two states (Berg 1995). 

From the experimentally obtained maps of network 
structure in the two states, a single network structure was 
derived that contained all vessels included in either of the 
observed maps. This network, which contains 125 vessels, 
is used as the basis for the present model. The observed 
maps included arterioles and capillaries up to the points 
at which they converged to form venules, but did not 
include information about the location of the venules. In 
order to complete the network structure, 12 venules were 
inserted, so as to connect all disconnected capillaries, 
while minimizing the length of venules added. This pro- 
cedure resulted in a 137-segment microvessel network. 
The pressure drop in the venules is typically a small frac- 
tion of the overall pressure drop in the network. More- 
over, the venules make a minor contribution to the 
overall oxygen exchange in the network. Therefore, 
the distributions of blood flow and oxygen transport in 
the network are relatively insensitive to the assumed 
geometry of this added venular network. 

Vessels in the 137-segment microvessel network are 
classified as arterioles, capillaries, or venules. Arterioles, 
which are capable of active regulation, are defined based 
on the branching pattern: If a vessel is the parent vessel 
at a diverging bifurcation and has a diameter in the 
dilated state of at least 8 fim, it is classified as an arteriole 
(Wiedeman et al. 1981). All other vessels (except the 
added venules) are classified as capillaries. This results in 




lOOmmHg 90 mmHg 66 mmHg 



15.4 mmHg 13.91 mmHg 12.91 mmHg 



Figure 1. Schematic of the network model. A complete flow pathway through the circulation is formed by an experimentally observed 
microvascular network, with the addition of upstream and downstream segments representing arteries, large arterioles, large venules, and 
veins. The distribution of intravascular pressures is shown for the reference state of the network model. 
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Figure 2. Network used in simulations. Arterioles (red) and 
capillaries (green) are derived from observations of hamster 
cremaster muscle (Berg 1995). Venules (blue) are added to form 
connections with downstream ends of capillaries. 



32 arterioles, 93 capillaries, and 12 venules (Fig. 2). Of 
the arterioles, six were terminal arterioles, giving about 15 
capillaries per terminal arteriole. This is similar to the 
number of capillaries previously stated to form a micro- 
vascular unit (Delashaw and Duling 1988). To form a 
complete flow pathway through the systemic circulation, 
two upstream vessels (artery A, and large arteriole LA) 
and four downstream vessels (large venules LV1, LV2, 
and veins VI, V2) are added, one of which (LA) is 
capable of active regulation, bringing the total number of 
vessels in the simulated network to 143. 



Network hemodynamics 

The microvessel network is represented as a set of 
interconnected segments, each with a defined diameter 
and length. Pressure-driven flow in each segment and 
conservation of mass at each node are assumed. Phase 
separation of RBCs at diverging bifurcations is imple- 
mented based on previously derived empirical relation- 
ships (Pries et al. 1996). Flow in each segment is 
governed by Poiseuille's law, with an apparent viscosity 



that depends on vessel diameter and includes effects of 
an endothelial surface layer (Pries and Secomb 2008). 
Using an iterative technique (Young 1954), the resulting 
nonlinear system of equations is solved at each time 
step with fixed diameters to obtain the apparent viscos- 
ity, flow rate, hematocrit, and wall shear stress in each 
segment. 

In the experimental observations (Berg 1995), vessels 
were considered to be flowing if they had observable RBC 
fluxes. As RBC flux was reported in steps of 25 cells/sec, 
the threshold for observable RBC flux is assumed to be 
^threshold = 12.5 cells/sec. In the following, vessels with 
RBC flux above this threshold are referred to as flowing 
and others as nonflowing. 

Flow regulation model 

The model for flow regulation is based on the modulation 
of arteriolar diameters according to the passive and active 
length-tension characteristics of VSM, and includes myo- 
genic, shear-dependent, and metabolic responses (Arciero 
et al. 2008; Carlson et al. 2008). Steady-state tension in 
the vessel wall (T total ) is represented as the sum of a pas- 
sive component and an active component generated by 
the VSM (Carlson and Secomb 2005): 



Ttotal — Tp ass + AT act 



(1) 



Here, A represents the local activation level of the mus- 
cle in each vessel, and has a range from 0 to 1, where 
A = 0 represents no vascular tone and A = 1 represents 
maximal vasoconstriction. The passive tension in the wall 
of an arteriole with diameter D is given by 

D 



T 

1 ry. 



Cpass exp 



c' 



Do 



(2) 



where C pass and C' pass are constants representing the magni- 
tude and diameter dependence of passive tension, and D 0 is 
the passive diameter at an intraluminal pressure of 
100 mmHg. The maximal active wall tension is described by 



C ac , exp 



C 



C" 

^act 



(3) 



where C act , C' act , and C" act are constants representing 
maximally active VSM peak tension, length dependence, 
and tension range, respectively. A target activation level is 
introduced, which is assumed to be a saturating function 
of the total vasoactive signal, S tone : 

Aotai = — \—z r • (4) 

1 + exp(-S tone ) 

The parameter S tone represents the combined effects of 
the myogenic, shear-dependent, and metabolic input 
stimuli (Arciero et al. 2008; Carlson et al. 2008): 
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•^tone C m y 0 T C snear T wa n 5 me ta ~T" C tone 7 (5) 

where C myo , C snea r> an d C" tone are constants representing 
the sensitivity of S tone to the various stimuli. The wall 
tension is given by T = PD/2, where P is the average seg- 
ment intraluminal pressure, T waU is the wall shear stress, 
and S meta is the conducted metabolic response signal, 
computed as described below. The model defined by 
equations (l)-(5) was shown to yield steady-state 
pressure-diameter relationships in close agreement with 
experimental observations of myogenic responses in 
arterioles of various tissues (Carlson and Secomb 2005; 
Carlson et al. 2008). 

Parameter values in the model for VSM tone were 
chosen to match those used previously (Berg 1995; Arciero 
et al. 2008; Carlson et al. 2008; Roy et al. 2012). For the 
parameter C act , the maximally active VSM peak tension, 
values were needed for smaller diameter vessels (<30 in 
the maximally dilated state) than those considered in the 
previous models. Therefore, C act was estimated by fitting a 
power law relationship to six data sets (Carlson and 
Secomb 2005) with the smallest average value of D u , 
including four data sets for hamster microvessels (Fig. 3). 
Values of parameters used in the model are given in 
Table 1. 

Oxygen transport 

Previous models for metabolic flow regulation considered 
oxygen exchange between each vessel and an associated 

2000 i 1 




D 0 (^m) 

Figure 3. Power law fit of C act versus D 0 for six experimental data 
sets (Carlson and Secomb 2005). Solid circles: experimental data. 
Curve: power law fit: C act = 1 .30 Do 1 ' 48 . In the fitting procedure, 
the data points were weighted by the number of vessels in each 
data set. 
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local tissue region (or "tissue sleeve") of fixed width sur- 
rounding the vessel (Arciero et al. 2008; Roy et al. 2012). 
A limitation of this type of oxygen transport model is 
that as blood flow rate approaches zero in an individual 
vessel, the partial pressure of oxygen (P0 2 ) must go to 
zero in that vessel. In reality, tissue cells may receive oxy- 
gen from any vessel that is sufficiently close. A decrease 
in flow in a particular vessel does not necessarily cause 
hypoxia if other well-oxygenated vessels are nearby. 
Therefore, models based on a fixed oxygen consumption 
rate per vessel length overestimate the flow-dependent 
changes in intravascular oxygen levels and exaggerate the 
ability of metabolic responses to maintain flow in all ves- 
sels. Such models are therefore not suitable for investigat- 
ing vascular recruitment. 

Here, a more realistic oxygen transport model is used. 
The model is capable of simulating a spatially varying 
oxygen field, including effects of diffusive interaction 
between all vessels in the network and the entire tissue 
region. The governing equation in the tissue is as follows: 

DaV 2 p = M, (6) 

where D and a are the diffusivity and solubility of oxygen 
and tissue, p is the P0 2 and M is the oxygen consump- 
tion rate, which is assumed to follow Michaelis-Menten 
kinetics: 



A Michaelis constant (P 0 ) of 1 mmHg is used (Secomb 
et al. 2004). Oxygen demand (M 0 ) is varied between 0.5 
and 8 cm 3 O 2 /100 cm 3 /min, corresponding to a range of 
metabolic conditions. Convective oxygen transport in ves- 
sels is simulated, including the effects of variable hemato- 
crit and the nonlinear oxygen binding characteristics of 
hemoglobin. Equation (6) is solved by a Green's function 
method (Secomb et al. 2004), in which each vessel is rep- 
resented as a set of discrete oxygen sources and the tissue 
is represented a set of discrete oxygen sinks. The oxygen 
field in the tissue is represented as a superposition of the 
fields resulting from the sources and sinks. This allows 
for more efficient computation, relative to other compu- 
tational approaches, as it reduces the number of 
unknowns needed to represent the oxygen field. 

The geometric locations of the added upstream and 
downstream vessels are not specified. For the purpose of 
modeling oxygen transport, it is assumed that the arterial 
vessels lying outside the region with observed network 
geometry (A and LA) are surrounded by a tissue sleeve of 
width 18.8 in which oxygen is consumed at a fixed 
rate. The sleeve width is chosen to correspond to a typical 
measured value of skeletal muscle capillary density of 
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Table 1. Parameter values used in the model. 



Description 


Parameter 


Value 


Unit 




Reference 


Oxygen demand 


Mo 


U.D LO O 


cm 3 O 2 /100 cm 3 /min 




Hill equation exponent 


n 


2. SB 






Klitzman et al. (1983) 


Half-maximal Hb saturation 


Pso 


26 


mmHg 




Klitzman et al. (1983) 


Time constant for diameter 


Td 


I 


sec 




Arciero et ai. \zuuo) 


Time constant for activation 




zu 


sec 




Koy ei ai. vzu i z; 


VSM activation tension sensitivity 


r- 

'-myo 


I . jjby/Uo 


cm/dyn 




Carlson and Secomb (2005) 


VSM activation shear stress sensitivity 


r- 

'-shear 


U.Uzdo 


cm 2 /dyn 




Lanson et ai. vzuuoj 


VSM activation conducted response 


Qneta 


5000 


\ i/jUM/cmj 




Koy et al. {Zu \ Z) 


sensitivity 












VSM constant 


(- tone 


— u.yuo to i Jz. i 






Carlson and Secomb (2005), 
Larison et al. (ZUUoj, 
Koy et ai. (zu I Z/ 


Passive tension strength 


r- 

*-pass 


D.DDD-iJo 


dyn/cm 




Carlson and Secomb (2005) 


Passive tension sensitivity 


pass 


-0.027-Do + 12.52 






Carlson et al. (2008) 


Maximally active VSM peak tension 


Cac, 


1.30D 0 148 


dyn/cm 




Carlson and Secomb (2005) 


Maximally active VSM length 


C'ac, 


-0.00059-Do + 0.773 






Carlson et al. (2008) 


dependence 












Maximally active VSM tension range 


C'act 


-0.00080-Do + 0.41 S 






Carlson et al. (2008) 


Passive reference vessel diameter 


Do 


8 to 20 






Berg (1995) 


Threshold for observable RBC flux 


F threshold 


12.5 


cells/sec 




Berg (1995) 


P0 2 where metabolic signal is 


Po 


1 


mmHg 




Secomb et al. (2004) 



half-maximal 



5007mm 2 (Arciero et al. 2008). Oxygen exchange by the 
added venous vessels (LV1, LV2, VI, and V2) is 
neglected. The Hill equation is used to calculate oxyhe- 
moglobin saturation as a function of vessel P0 2 , with 
parameters of n = 2.55 and P 50 = 26 mmHg, based on 
estimates in hamster cremaster muscle (Klitzman et al. 
1983). Oxygen content of all inflowing vessels to the net- 
works is prescribed in terms of the P0 2 - At the arterial 
inflow to the network, a P0 2 of 100 mmHg is assumed. 
The model includes oxygen loss from arterioles (Duling 
and Berne 1970). At venular inflows to the network, a 
P0 2 of 20 mmHg is assumed based on estimates of venu- 
lar saturation (Ellsworth et al. 1988). 

Conducted metabolic signal 

The origins of the metabolic signals involved in local flow 
regulation are not definitely established. Roy et al. (2012) 
analyzed flow regulation in heterogeneous microvascular 
networks, assuming either wall-derived or RBC-derived 
metabolic signals, and showed that a wall-derived signal 
was more effective in matching perfusion to local oxygen 
demand. In reality, several mechanisms likely contribute 
to metabolic flow regulation, involving wall-derived, tis- 
sue-derived, and/or RBC-derived signals. The involvement 
of mechanisms independent of RBCs has been shown 
(Ngo et al. 2010). In this study, the metabolic signals in 
each segment are assumed to be generated by signals from 



the vessel walls, as this type of signaling was found to be 
effective for local flow regulation in heterogeneous net- 
works (Roy et al. 2012). Although several potential mech- 
anisms leading to the generation of wall-derived signals 
are known, their characteristics remain incompletely 
understood, including their ranges of P0 2 sensitivity and 
their individual contributions to metabolic regulation 
under various conditions. Here, it is assumed that the 
local signal S; oc reflects the local oxygen deficit and is pro- 
portional to the difference between oxygen demand (M 0 ) 
and oxygen consumption (M, given by eq. 7), that is, 



where C meta is a parameter describing the strength of the 
signal and P 0 = 1 mmHg, corresponding to the Michaelis 
constant for tissue oxygen consumption. 

The signal is assumed to be conducted upstream to the 
arterioles, resulting in metabolic control of arteriolar tone 
(Segal et al. 1989), with exponential decay according to a 
characteristic length constant. Effects of varying the length 
constant were tested in the model and shown to have lit- 
tle effect, as long as it was not greatly decreased. At bifur- 
cations diverging in the upstream direction, the 
conducted signal is partitioned in proportion to the vessel 
diameters; at upstream converging bifurcations, the signal 
is summed (Roy et al. 2012). The total metabolic signal 
in each vessel, 

Smeta> is the sum of the propagated signal 
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from downstream vessels and the signal generated locally 
(Sjoc). Virtual segments are added distal to draining ves- 
sels to generate conducted signals at the network out- 
flows. 



Dynamics of arterioles 

The dynamic behavior of each arteriole is described by 
differential equations for the variation with time of diam- 
eter and activation (Arciero et al. 2008; Carlson et al. 
2008; Roy et al. 2012): 



dD \D C , 

— = (T - 

dt T rf T c 



dA 
It 



1 



* (Atotal ~ A), 



(9) 
(10) 



where D c and T c are the diameter and wall tension in the 
reference state (described below), and r d and r a are time 
constants. In these simulations, the eventual behavior of 
the system is analyzed, not the transient response, so the 
assumed values of the time constants do not generally affect 
the results. Available data indicate that T d is around 1 sec, 
and r a is between 1 sec and 1 min, so an intermediate value 
(20 sec) was used (Roy et al. 2012). Equations (9) and (10) 
are integrated from t = 0 to 200 sec using an explicit Euler 
method. After 200 sec, the system is found either to reach a 
steady state or to show stable oscillations, consistent with 
occurrence of vasomotion (Arciero and Secomb 2012). For 
purposes of analysis, the final values of system variables are 
defined by averaging over the interval from 100 to 200 sec. 
The simulation of oxygen transport is the most computa- 
tionally demanding part of the calculation. The iterative 
method used in the Green's function method at each time 
step involves solutions of large linear systems, which can be 
implemented with parallel algorithms. Simulation of 
200 sec of real time requires 800 computational time steps, 
with each time step taking 10 to 60 sec of computer time 
using a graphical processing unit-based parallel processing 
system. 

Reference state 

A reference state is defined for the network with a moder- 
ate level of arteriolar tone, corresponding to a relatively 
low level of oxygen consumption in the skeletal muscle. 
The distribution of pressure drops in the network model 
is chosen to correspond where possible to the previous 
model (Arciero et al. 2008), with pressure drops in the A, 
LV, and V segments of 10, 1.49, and 1 mmHg. The pres- 
sure drop across the microvessel network is chosen to 
minimize the total squared variation between vessel RBC 
fluxes in the reference state and those observed experi- 



mentally data in the control state (Berg 1995), resulting 
in a pressure of 66 mmHg at the arterial side of the mi- 
crovessel network. This was achieved by adjusting the 
length of the LA segment. The arterial inflow and venous 
outflow pressures are set to 100 and 12.91 mmHg, respec- 
tively (Roy et al. 2012). The resulting distribution of pres- 
sures is indicated in Figure 1. 

To establish the distribution of tone in the reference 
state, values of S tone in each arteriole are determined as 
follows. If an arteriole had observable RBC flux in the 
maximally dilated experimental network (Berg 1995), 
then the passive vessel diameter at a pressure of 
100 mmHg (D 0 ) is set as the measured diameter in the 
maximally dilated state. If an arteriole had observable 
RBC flux also in the control experimental network (Berg 
1995), then the diameter in the model reference state 
(D c ) is set as the measured diameter in the control state. 
The reference state is assumed to be at equilibrium, so 
that T = T tota i and A = A tota i, which allows for an explicit 
calculation of S tone in each arteriole based on equa- 
tions (1) and (4): 



Stone — l n 



pn _ 2T 



9 ( 7 1 _1_ Tmax A 



PD, 



(11) 



For some vessels in the experimental data set, the 
diameter in the maximally dilated state was less than or 
equal to the diameter in the control state. Such behavior 
might result from a decrease in pressure in some arteri- 
oles due to overall vasodilation. In these cases, it is 
assumed that the vessels are close to maximal dilation in 
the reference state; that is, D c w D 0 , and S tone is set to a 
large negative value (—1000) so that A « 0. Some arteri- 
oles had observable RBC fluxes only in the dilated state 
of the network. In the reference state, these arterioles are 
assigned small diameters D c in the range 2.9-3.8 /im, such 
that all have RBC fluxes below ^threshold- 

Of the capillaries, 32 were observed as flowing only in 
the dilated state, implying that they should have RBC 
fluxes less than ^threshold m the reference state, but greater 
than threshold when arteriolar diameters are set to their 
dilated values. When arteriolar diameters were initially set 
to their D c values as described above, only 11 of the 32 
capillaries had RBC fluxes less than threshold in the simu- 
lated reference state. The absolute root mean square 
(RMS) error in individual diameter measurements in 
microvessel networks was estimated as ±1.2 jim based on 
imaging limitations (Pries et al. 1994). Therefore, small 
adjustments (<1.2-^m increase or decrease) in capillary 
diameter were made where this led to an increase in the 
number of nonflowing capillaries in the reference state. 
After these changes were made, 19 of the 32 capillaries 
flowing only in the dilated state are nonflowing in the 
reference state, and all 32 are flowing when arteriolar 
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diameters are set to their dilated values. In total, of the 
125 capillaries and arterioles in the microvessel network, 
27 (21.6%) are nonflowing in the reference state, and 
therefore have the capacity for recruitment. 

As is evident from the above description, the procedure 
for defining the control state involves a number of 
assumptions. The resulting reference state may not accu- 
rately represent, at the level of individual vessels, the con- 
ditions that existed in the experimental preparation. 
However, the procedure used was designed to ensure that 
distributions of geometric and hemodynamic parameters 
and their degree of heterogeneity in the reference configu- 
ration are similar to those existing in the in vivo prepara- 
tion. The purpose of the reference state is to provide an 
internally consistent network with realistic properties as a 
reference configuration, which serves as a control for fur- 
ther simulations in which hemodynamic or metabolic 
parameters are varied. 



120 



100 



80- 



60 



40 



Arterioles 
Capillaries 




Experimental observations 



Model: flowing vessels 



Figure 4. Numbers of vessels in experimentally observed networks 
and numbers of flowing vessels in model simulations. 



Results 

The numbers of flowing arterioles and capillaries observed 
experimentally and predicted by the model are shown in 
Figure 4. In the experimental observations, the total num- 
ber of flowing vessels increased by 26% (from 85 to 107) 
between the control and vasodilated states. The number of 
flowing arterioles and capillaries in the dilated state (107) is 
fewer than the total number of such vessels in the com- 
bined network consisting of the union of the networks in 
the control and vasodilated states (125) because a number 
of vessels (18) were observed to be flowing in the control 
state but not in the dilated state. Such behavior may result 
from hemodynamic effects including redistribution of 
hematocrit with vasodilation, or it may reflect incomplete 
visualization of all flowing vessels in the network in the 
experiments. The combined network is used as the basis for 
the model simulations. 

As shown in Figure 4, the model predicts that the 
numbers of flowing arterioles and capillaries increase 
progressively with increasing oxygen demand. This sup- 
ports the hypothesis that vessel recruitment can occur as 
a consequence of local metabolic regulation of blood 
flow in a network with heterogeneous structure. Overall, 
the number of flowing arterioles and capillaries increases 
by about 16% (from 98 to 114), as demand (M 0 ) is 
raised from 1 to 8 cm 3 O 2 /100 cm 3 /min. The spatial dis- 
tributions of flowing and nonflowing vessels are shown 
in Figure 5 for oxygen demand M 0 of 1, 2, 4, and 8 cm 3 
O 2 /100 cm 3 /min. In general, increasing oxygen demand 
causes transitions from nonflowing to flowing states, but 
a few vessels stop flow with increasing demand (Fig. 5A). 
In many instances, recruitment or derecruitment of 
individual capillaries, rather than groups fed by a single 



terminal arteriole, is predicted with changes in oxygen 
demand. 

The effects on predicted tissue oxygenation of including 
flow regulation and hematocrit partition in the model are 
illustrated in Figure 6, as spatial distributions (upper pan- 
els) and frequency distributions (lower panels). For the 
reference state, with oxygen demand M 0 = 1 cm 3 0 2 / 
100 cm 3 /min, virtually no hypoxic tissue is predicted 
(Fig. 6A). When oxygen demand is increased to M 0 = 4 
cm 3 O 2 /100 cm 3 /min, without metabolic flow regulation 
(i.e., S meta = 0 in each vessel) or phase separation 
included in the model, the distribution of tissue P0 2 is 
shifted toward lower levels (Fig. 6B). When metabolic 
flow regulation is included (Fig. 6C), the tissue P0 2 dis- 
tribution shifts toward higher levels compared to the case 
without metabolic flow regulation. Figure 6D shows the 
distribution with both metabolic flow regulation and 
phase separation included. The inclusion of nonuniform 
hematocrit distribution further shifts the tissue P0 2 dis- 
tribution to higher levels. 

In Figure 7, variables reflecting blood flow and tissue 
oxygenation are plotted as functions of oxygen demand 
from M 0 = 0.5 to 8 cm 3 O 2 /100 cnrVmin for three cases: 
constant hematocrit without metabolic regulation, con- 
stant hematocrit with metabolic regulation, and variable 
hematocrit with metabolic regulation. Figure 7A shows 
total blood flow to the network. With metabolic regula- 
tion, inflow rate increases as oxygen demand is increased. 
Introducing phase separation has little effect on total 
flow. Figure 7B shows the number of flowing arterioles 
and capillaries. The number of vessels with RBC flux 
above threshold is lower when phase separation is 
included. The number increases 13% between the lowest 
and highest oxygen demand states with constant hemato- 
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Figure 5. Spatial distributions of nonflowing vessels and flowing arterioles, capillaries, and venules, when (A) M 0 = 1; (B) M 0 = 2; (C) M 0 = 4; 
(D) M a = 8 cm 3 O 2 /100 cm 3 /min. Black arrows indicate vessels that start flow with an increase in oxygen demand. Black arrows with an 
asterisk in panel A indicate vessels that stop flow with an increase in oxygen demand. 



exit, whereas the increase is 25% with phase separation. 
Thus, phase separation increases the extent of vessel 
recruitment as oxygen demand is increased. As expected, 
the network inflow rate and the number of flowing vessels 
are independent of oxygen demand when metabolic flow 
regulation is not included (Fig. 7 A and B). 

Median tissue P0 2 declines with increasing oxygen 
demand (Fig. 7D). Inclusion of metabolic flow regulation 
leads to an increase of about 10 mmHg in median tissue 
P0 2 - The increase is larger when phase separation is 
included, relative to the uniform hematocrit case. These 
differences are present for the entire range of oxygen 
demand considered. Figure 7C shows the hypoxic fraction 
(fraction of tissue with P0 2 < 1 mmHg). For oxygen 
demand up to M 0 = 2 cm 3 O 2 /100 cm 3 /min, there is less 
than 1% hypoxia in all cases. However, for M 0 = 4, with 
no metabolic flow regulation, the hypoxic fractions are 
0.031 and 0.028 for variable H D (not shown) and con- 
stant H D , respectively. When metabolic flow regulation is 
included, the hypoxic fraction is much lower, with values 
of 0.006 and 0.005, respectively. Similarly, for M 0 = 8, the 
hypoxic fractions without metabolic flow regulation 
increase to 0.14 (for both variable and constant H D ) com- 
pared to 0.06 (variable H D ) and 0.061 (constant H D ) with 
metabolic regulation included. Thus, for M 0 > 2, tissue 
hypoxic fraction is much lower when metabolic flow reg- 
ulation is included, as would be expected. 

Discussion 

Metabolic regulation is necessary for matching blood flow 
to widely varying oxygen demand in striated muscle. 



Capillary recruitment, defined as variation in the number 
of capillaries containing flowing RBCs, is observed to 
occur with changes in metabolic demands and is consid- 
ered to contribute to metabolic control of blood flow. 
This study uses a theoretical model to test the hypothesis 
that vessel recruitment can occur as a consequence of 
local changes in arteriolar tone, acting in combination 
with phase separation at diverging microvessel bifurca- 
tions. The results support the hypothesis. When applied 
to a heterogeneous network structure in striated muscle 
derived from a hamster cremaster preparation (Berg 
1995), the model predicts that metabolic vasodilation not 
only increases overall flow rate (Fig. 7 A) but also 
increases the number of flowing vessels (Figs. 4, 5, and 
7B) defined in terms of RBC flux above a threshold value. 
With increasing oxygen demand, vessel recruitment is 
predicted at both the capillary and arteriolar levels. 

As discussed in the Introduction, the mechanism of 
recruitment has been unclear. The most distal site for 
control of capillary perfusion in skeletal muscle is the ter- 
minal arteriole, which controls flow to a microvascular 
unit comprising multiple capillaries (around 20), whereas 
recruitment of individual capillaries is observed. This 
study provides a potential resolution to these apparently 
conflicting observations. In the simulations, recruitment 
and derecruitment of individual capillaries are predicted 
to occur as a consequence of arteriolar dilation and con- 
striction. Here, as in experimental studies, capillary perfu- 
sion is defined in terms of RBC flux above a minimum 
observable threshold. In a heterogeneous network struc- 
ture, the unequal partition of hematocrit at diverging 
bifurcations causes vessel hematocrits to vary depending 
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Figure 6. Contour plots (top panels) and histograms (bottom panels) of tissue P0 2 distribution at hundreds of tissue points. Oxygen levels in 
the vessels and in the surrounding tissue are color-coded according to the scale at left (in mmHg). (A) M 0 = 1 cm 3 O 2 /100 cm 3 /min. 
(B) M 0 = 4 cm 3 O 2 /100 cm 3 /min, with constant hematocrit and without metabolic flow regulation. (C) M 0 = 4 cm 3 O2/IOO cm 3 /min, with 
constant hematocrit and with metabolic flow regulation. (D) M 0 = 4 cm 3 O 2 /100 cm 3 /min, with nonuniform hematocrit and with metabolic 
flow regulation. 



on the flow distribution. As a consequence, recruitment 
of capillaries within a single microvascular unit can occur 
sequentially, rather than synchronously, during gradual 
vasodilation of the feeding arteriole. According to this 
concept, modulation of flow and vascular recruitment are 
two aspects of the same physiological response, and it is 
not meaningful to quantify their individual contributions 
to the metabolic regulation of flow. 

The effects of metabolic flow regulation, with and 
without unequal hematocrit partition at vessel bifurca- 
tions, on predicted hemodynamic and oxygen transport 
characteristics are illustrated in Figure 7. Including phase 
separation has little effect on the total flow to the net- 
work (Fig. 7A), as variations in individual vessel flow 
rates tend to cancel when aggregated. However, phase 
separation has a substantial effect on the number of flow- 
ing vessels (Fig. 7B). With phase separation, a larger 
number of vessels have low hematocrit and low RBC 
fluxes, and are classified as nonflowing. The extent of 



recruitment and derecruitment is greater with phase sep- 
aration because decreases in either flow rate or vessel 
hematocrit may cause the RBC flux to decrease below 
the threshold level. 

These dual mechanisms contributing to the distribution 
of RBC flux in microvascular networks and the recruit- 
ment and derecruitment of capillaries are illustrated sche- 
matically in Figure 8. Regardless of whether phase 
separation is included in the model, low arteriolar flows 
may result in capillaries with RBC fluxes below ^threshold 
that can be recruited when arteriolar flow increases in 
response to increased oxygen demand. A second mecha- 
nism becomes evident itself only when phase separation is 
included. Capillary RBC fluxes below -F t h r eshoid ma Y occur 
because a parent arteriole receives a low hematocrit. With 
vasodilation, the arteriole may receive a higher hemato- 
crit, leading to recruitment of the capillaries that it feeds. 
Capillary recruitment is thereby augmented as a result of 
unequal hematocrit partition at vessel bifurcations. 
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Figure 7. Variables describing network blood flow and tissue oxygenation, as functions of oxygen demand. (A) Relative arteriolar inflow rate. 
(B) Relative number of flowing vessels. (C) Hypoxic tissue fraction (defined as P0 2 < 1 mmHg). (D) Median tissue P0 2 . (•) Variable H D (with 
phase separation) with metabolic flow regulation, (o) Constant H D (no phase separation) with metabolic flow regulation. (T) Constant H D 
without metabolic flow regulation. 



An unexpected model prediction is that inclusion of 
phase separation in the model with metabolic flow reg- 
ulation leads to improved tissue oxygenation, relative 
to the behavior without phase separation (i.e., with 
constant hematocrit), as shown by the shift to higher 
values in the distribution of P0 2 (Fig. 6C and D) and 
the median tissue P0 2 (Fig. 7D). This improvement in 
oxygenation occurs despite the fact that fewer vessels 
are flowing when phase separation is included 
(Fig. 7B). Phase separation is thus seen to contribute 
to metabolic flow regulation. By directing increased 
hematocrit preferentially to vessels with increased flow 
rates, it provides a mechanism to amplify the changes 
in RBC flux resulting from modulation of arteriolar 
diameters. In the present simulations, the inclusion of 
phase separation has only slight effects on the fraction 
of hypoxic tissue (Fig. 7C). However, the hypoxic 
regions are small and lie mainly near the boundaries of 



the network (see Fig. 6), and so this result may not be 
representative of the behavior that would occur with 
more widespread hypoxia. 

In the present model, the signal for metabolic regula- 
tion in response to changes in oxygen levels is assumed to 
be generated in the vessel wall, according to equation (8). 
A previous theoretical study (Roy et al. 2012) examined 
the distribution of blood flow and tissue oxygenation in a 
heterogeneous network, assuming that the signal for met- 
abolic flow regulation originated within RBCs. In that 
case, inclusion in the model of phase separation resulted 
in poorer oxygenation because a reduction in RBC flux to 
a given vessel generally caused a reduction in the RBC- 
derived metabolic signal, thereby exacerbating the reduc- 
tion in flow. Taken together, these models show that 
phase separation can have a substantial effect on local 
metabolic flow regulation, but its effect depends on the 
metabolic signaling mechanisms involved. 
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Figure 8. Schematic diagram showing two mechanisms for 
capillary recruitment and derecruitment in response to changes in 
arteriolar diameter. Black arrows indicate magnitude of RBC flux. 
Gray arrows indicate arteriolar vasoconstriction. Gray areas within 
vessels represent regions containing RBCs. (A) Initial configuration 
with dilated arteriole. (B) Effect of vasoconstriction without phase 
separation. Derecruitment may occur because RBC flux drops below 
threshold level. (C) Effect of vasoconstriction with phase separation. 
Derecruitment may occur because hematocrit entering the arteriole 
drops to a low level or zero. 



Estimates of the capacity for capillary recruitment in 
skeletal muscle vary widely (Clark et al. 2008). In some 
experiments, it has been found to be at least 100% and 
possibly as high as 200%, relative to the number of flow- 
ing capillaries at rest (Honig et al. 1982; Lindbom 1983). 
In the present simulations, almost 80% of the arterioles 
and capillaries in the network are flowing in the reference 
state, and so the capacity for recruitment is limited to 
about 25%. Furthermore, the range of flow regulation in 
this model study is restricted according to the observed 
changes in diameter between the two states of the cremas- 
ter preparation (Berg 1995). The physiological range of 
flow rates in the cremaster muscle is not as wide as that 
in skeletal muscle (Segal 2005). The geometry of the cre- 
master muscle, with layers of fibers running in different 
directions, differs from that of skeletal muscle. Despite 
these differences, there is no reason to expect that the 
recruitment mechanism described here would not also be 
operative in skeletal muscle. The ranges of flow modula- 
tion and capillary recruitment possible with the hypothe- 
sized mechanisms have no inherent restrictions. Given a 
network with more vasoconstriction and fewer flowing 
capillaries in the resting state, the present model would 
predict a wider range of flow regulation and vascular 
recruitment. 
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